Analysis of Hurricanes in the Atlantic Basin

Background

The HURDAT dataset is provided by the NOAA, containing data in six-hour time intervals for all hurricanes and tropical storms since 1851. For the purposes of this analysis, we will look at all hurricanes in the North Atlantic Ocean, affecting North America and the Caribbean.

The data is provided as a text file, containing a multidimensional array for storms and individual datapoints for each time interval. The purpose of this notebook is to parse this raw data into cleanly formatted and labeled data that can be used for data exploration and analysis.

This data cleaning will include creating features for the tropical storm/hurricane as a whole (e.g. did it make landfall/what was the max category/duration of system).

Contents

1. Parse Hurdat Dataset
   1a. Load Data
   1b. Format Into Dictionary
   1c. Feature Engineering for Individual Storms
   

2. Analyze Hurricane Paths

3. Exploratory Data Analysis of Atlantic Hurricanes 1850-2017
   3a. Number of Storms
   3b. Storm Duration Over Time

1. Parse Hurdat Dataset

Data Description: https://www.nhc.noaa.gov/data/hurdat/hurdat2-format-atlantic.pdf

import os
import re
import pandas as pd
import numpy as np
import datetime
import json
import folium
import seaborn as sns
import chart_studio.plotly as py
import plotly.graph_objs as go
from urllib.request import urlretrieve
from tqdm import tqdm_notebook as tqdm
from plotly.tools import FigureFactory as FF
from plotly import offline
from plotly.offline import iplot
from plotly.offline import init_notebook_mode, iplot
init_notebook_mode()

1a. Load Data

Hurricane data from the NOAA comes in a pretty raw format. We will first:

1. Read in the raw data using read_hurdat()
2. Convert the lat/lon which is in a string format using directional indivators into numeric values using convert_lat_lon()
3. Categorize the wind speed severity of a hurricane at any given moment using hurr_category()

def read_hurdat(url, local_fname, location):
    """
    Loads and parses the HURDAT dataset from the NOAA website, given
    the input url and location to store the parsed file locally.
    """
    if not os.path.exists(local_fname):
        urlretrieve(url, local_fname)

    records = []
    with open(local_fname,'r') as f:
        for line in f:
            if line.startswith(location):
                record = line.strip()
                reports = []
                records.append((record, reports))
            else:
                reports.append(line.strip())
                
    return records

def convert_lat_lon(value, col):
    """
    Lat/lon is encoded with the numeric value plus E/W/N/S. We want to 
    convert this into an absolute decimal value for maps.
    """
    if col=='lon':
        amount = -float(re.sub('[EW]', '', value)) if 'W' in value else float(re.sub('[EW]', '', value))
    elif col=='lat':
        amount = -float(re.sub('[NS]', '', value)) if 'S' in value else float(re.sub('[NS]', '', value))
    
    return amount

def hurr_category(max_wind):
    """
    Defining storm category based on maximum observed winspeed at the time.
    """
    max_wind = int(max_wind)
    if max_wind <= 73:
        cat = 'TS'
    if (74 <= max_wind <= 95):
        cat = 1
    if (96 <= max_wind <= 110):
        cat = 2
    if (111 <= max_wind <= 129):
        cat = 3
    if (130 <= max_wind <= 156):
        cat = 4
    if max_wind >= 157:
        cat = 5
    return cat

url = "https://www.nhc.noaa.gov/data/hurdat/hurdat2-1851-2017-050118.txt"
local_fname = "../data/hurdat2.txt"

records = read_hurdat(url, local_fname, "AL") # AL for atlantic hurricanes

Lets look at what a raw record looks like (Hurricane RINA):

records[-5]

('AL152017,              MARIA,     68,',
 ['20170916, 1200,  , TD, 12.2N,  49.7W,  30, 1006,    0,    0,    0,    0,    0,    0,    0,    0,    0,    0,    0,    0,',
  '20170916, 1800,  , TS, 12.2N,  51.7W,  40, 1004,   40,    0,    0,   40,    0,    0,    0,    0,    0,    0,    0,    0,',
  '20170917, 0000,  , TS, 12.4N,  53.1W,  45, 1002,   40,   30,    0,   40,    0,    0,    0,    0,    0,    0,    0,    0,',
  '20170917, 0600,  , TS, 12.8N,  54.4W,  55,  994,   50,   40,    0,   50,   20,   20,    0,   20,    0,    0,    0,    0,',
  '20170917, 1200,  , TS, 13.3N,  55.7W,  60,  990,   60,   40,   30,   50,   30,   20,    0,   20,    0,    0,    0,    0,',
  '20170917, 1800,  , HU, 13.6N,  57.0W,  65,  986,   70,   60,   40,   60,   30,   20,    0,   20,   15,    0,    0,   10,',
  '20170918, 0000,  , HU, 14.0N,  58.0W,  75,  979,   90,   60,   40,   70,   30,   20,   20,   20,   15,   10,   10,   10,',
  '20170918, 0600,  , HU, 14.3N,  59.0W,  80,  977,   90,   60,   40,   70,   30,   20,   20,   30,   15,   10,   10,   10,',
  '20170918, 1200,  , HU, 14.5N,  59.7W, 100,  967,  110,   90,   70,   90,   30,   30,   30,   30,   15,   15,   15,   15,',
  '20170918, 1800,  , HU, 14.9N,  60.4W, 110,  956,  110,   90,   80,   90,   30,   30,   30,   30,   20,   15,   15,   20,',
  '20170919, 0000,  , HU, 15.3N,  61.1W, 145,  924,  110,  110,   80,   90,   40,   40,   30,   40,   25,   25,   20,   25,',
  '20170919, 0115, L, HU, 15.4N,  61.3W, 145,  922,  110,  110,   80,   90,   40,   40,   30,   40,   25,   25,   20,   25,',
  '20170919, 0600,  , HU, 15.7N,  61.9W, 135,  940,  110,  110,   80,  100,   50,   50,   30,   50,   25,   25,   20,   25,',
  '20170919, 1200,  , HU, 16.1N,  62.7W, 140,  931,  120,  110,   80,  100,   80,   50,   40,   80,   25,   25,   20,   30,',
  '20170919, 1800,  , HU, 16.6N,  63.5W, 145,  920,  120,  110,   80,  100,   80,   60,   60,   80,   30,   30,   30,   30,',
  '20170920, 0000,  , HU, 17.0N,  64.3W, 150,  909,  130,  110,  100,  110,   80,   70,   60,   80,   50,   45,   35,   40,',
  '20170920, 0300, I, HU, 17.3N,  64.7W, 150,  908,  130,  110,  100,  110,   80,   70,   60,   80,   50,   45,   35,   40,',
  '20170920, 0600,  , HU, 17.6N,  65.1W, 140,  913,  130,  110,  100,  110,   80,   70,   60,   80,   50,   45,   35,   40,',
  '20170920, 1015, L, HU, 18.0N,  65.9W, 135,  920,  130,  110,  100,  110,   80,   70,   60,   80,   50,   45,   35,   40,',
  '20170920, 1200,  , HU, 18.2N,  66.2W, 115,  935,  130,  110,  100,  110,   80,   70,   60,   80,   50,   45,   35,   40,',
  '20170920, 1800,  , HU, 18.6N,  67.0W,  95,  959,  130,  110,  100,  110,   70,   70,   60,   70,   50,   45,   35,   40,',
  '20170921, 0000,  , HU, 19.0N,  67.6W,  95,  958,  130,  110,  110,  110,   70,   70,   60,   70,   50,   30,   40,   40,',
  '20170921, 0600,  , HU, 19.4N,  68.2W, 100,  959,  130,  110,  110,  110,   70,   70,   60,   70,   50,   30,   40,   40,',
  '20170921, 1200,  , HU, 19.9N,  68.8W, 100,  959,  130,  110,  110,  110,   80,   80,   60,   80,   50,   50,   40,   50,',
  '20170921, 1800,  , HU, 20.5N,  69.5W, 105,  960,  140,  130,  110,  110,   90,   80,   60,   80,   60,   50,   40,   50,',
  '20170922, 0000,  , HU, 20.8N,  70.0W, 110,  953,  140,  120,  100,  120,   90,   80,   60,   80,   60,   50,   40,   50,',
  '20170922, 0600,  , HU, 21.2N,  70.5W, 110,  959,  140,  120,  100,  130,   90,   80,   60,   80,   60,   50,   40,   50,',
  '20170922, 1200,  , HU, 21.9N,  70.9W, 110,  958,  140,  120,  100,  130,   90,   80,   60,   90,   60,   40,   40,   40,',
  '20170922, 1800,  , HU, 22.8N,  71.2W, 110,  959,  140,  130,  100,  120,   80,   60,   60,   70,   50,   40,   35,   30,',
  '20170923, 0000,  , HU, 23.7N,  71.6W, 105,  953,  160,  130,  100,  120,   80,   60,   60,   70,   50,   40,   35,   30,',
  '20170923, 0600,  , HU, 24.4N,  71.9W, 100,  952,  170,  170,  100,  120,   90,   90,   60,   70,   50,   40,   35,   30,',
  '20170923, 1200,  , HU, 25.1N,  72.1W, 100,  952,  170,  170,  100,  150,   90,   90,   60,   70,   50,   40,   35,   40,',
  '20170923, 1800,  , HU, 25.9N,  72.3W, 100,  952,  210,  200,  100,  150,   90,   80,   60,   80,   50,   40,   35,   45,',
  '20170924, 0000,  , HU, 26.6N,  72.4W, 100,  945,  210,  200,  100,  150,  100,   80,   60,   80,   50,   40,   35,   45,',
  '20170924, 0600,  , HU, 27.5N,  72.6W,  95,  942,  210,  200,  120,  160,  100,   80,   60,   80,   50,   40,   35,   45,',
  '20170924, 1200,  , HU, 28.4N,  72.8W,  95,  947,  210,  200,  120,  160,  100,   80,   60,   80,   50,   40,   35,   45,',
  '20170924, 1800,  , HU, 29.1N,  72.9W,  90,  943,  200,  200,  140,  160,  100,  100,   70,   90,   50,   40,   35,   45,',
  '20170925, 0000,  , HU, 29.7N,  72.9W,  85,  947,  200,  200,  140,  160,  100,  100,   70,   90,   50,   40,   35,   45,',
  '20170925, 0600,  , HU, 30.3N,  72.9W,  75,  954,  200,  200,  160,  170,  100,  100,   90,   90,   60,   60,    0,    0,',
  '20170925, 1200,  , HU, 30.8N,  73.0W,  70,  961,  200,  200,  160,  160,  120,  120,   90,   90,   80,   80,    0,    0,',
  '20170925, 1800,  , HU, 31.4N,  73.1W,  70,  966,  180,  180,  160,  160,  120,  120,   90,   90,   90,   90,    0,    0,',
  '20170926, 0000,  , HU, 32.0N,  73.1W,  70,  966,  180,  180,  160,  160,  120,  120,   90,   90,   90,   90,    0,    0,',
  '20170926, 0600,  , HU, 32.6N,  73.1W,  65,  970,  200,  210,  160,  160,  120,  120,   90,   90,   90,   90,    0,    0,',
  '20170926, 1200,  , HU, 33.3N,  73.1W,  65,  970,  200,  210,  160,  160,  120,  120,   90,   90,   90,   90,    0,    0,',
  '20170926, 1800,  , HU, 33.9N,  73.1W,  65,  975,  200,  200,  160,  160,  120,  120,   80,  100,   90,   90,    0,   90,',
  '20170927, 0000,  , HU, 34.4N,  73.0W,  65,  975,  200,  200,  160,  160,  120,  120,   80,  100,   90,   90,    0,   90,',
  '20170927, 0600,  , HU, 34.9N,  72.9W,  65,  976,  200,  200,  160,  150,  120,  120,    0,  100,   90,    0,    0,   90,',
  '20170927, 1200,  , HU, 35.4N,  72.8W,  65,  977,  200,  200,  160,  150,  120,  120,    0,  100,   90,    0,    0,   90,',
  '20170927, 1800,  , HU, 36.0N,  72.6W,  65,  979,  200,  200,  160,  150,  120,  120,    0,  100,   90,    0,    0,   90,',
  '20170928, 0000,  , HU, 36.6N,  72.2W,  65,  979,  200,  200,  160,  150,  120,  120,    0,  100,   90,    0,    0,   90,',
  '20170928, 0600,  , TS, 36.7N,  71.3W,  60,  982,  180,  210,  160,  150,   90,   90,   70,   80,    0,    0,    0,    0,',
  '20170928, 1200,  , TS, 36.8N,  70.0W,  60,  982,  180,  210,  160,  150,   90,   90,   70,   80,    0,    0,    0,    0,',
  '20170928, 1800,  , TS, 36.8N,  68.6W,  55,  985,  180,  210,  160,  150,   90,   90,   70,   80,    0,    0,    0,    0,',
  '20170929, 0000,  , TS, 36.9N,  66.8W,  55,  985,  180,  200,  160,  150,   90,   90,   70,   80,    0,    0,    0,    0,',
  '20170929, 0600,  , TS, 37.0N,  64.6W,  50,  987,  100,  210,  220,  180,    0,   90,  100,    0,    0,    0,    0,    0,',
  '20170929, 1200,  , TS, 37.0N,  62.0W,  50,  988,  100,  210,  220,  180,    0,   90,  100,    0,    0,    0,    0,    0,',
  '20170929, 1800,  , TS, 37.4N,  59.0W,  50,  988,  100,  220,  220,  180,    0,  100,  100,    0,    0,    0,    0,    0,',
  '20170930, 0000,  , TS, 38.1N,  55.6W,  50,  988,  100,  220,  220,  180,    0,  100,  100,    0,    0,    0,    0,    0,',
  '20170930, 0600,  , TS, 39.1N,  52.2W,  50,  988,  110,  200,  200,  180,    0,  100,  100,    0,    0,    0,    0,    0,',
  '20170930, 1200,  , TS, 40.0N,  48.8W,  50,  988,  110,  200,  200,  180,    0,  100,  100,    0,    0,    0,    0,    0,',
  '20170930, 1800,  , EX, 41.2N,  45.6W,  45,  991,  110,  200,  200,  150,    0,    0,    0,    0,    0,    0,    0,    0,',
  '20171001, 0000,  , EX, 42.2N,  42.6W,  45,  994,  110,  200,  200,  150,    0,    0,    0,    0,    0,    0,    0,    0,',
  '20171001, 0600,  , EX, 43.4N,  39.4W,  45,  996,  110,  200,  200,  150,    0,    0,    0,    0,    0,    0,    0,    0,',
  '20171001, 1200,  , EX, 44.9N,  35.5W,  45,  999,  110,  200,  200,  150,    0,    0,    0,    0,    0,    0,    0,    0,',
  '20171001, 1800,  , EX, 46.5N,  31.0W,  45, 1003,  110,  200,  200,  150,    0,    0,    0,    0,    0,    0,    0,    0,',
  '20171002, 0000,  , EX, 47.5N,  26.5W,  40, 1005,    0,    0,  200,  150,    0,    0,    0,    0,    0,    0,    0,    0,',
  '20171002, 0600,  , EX, 48.0N,  22.0W,  40, 1012,    0,    0,  200,  150,    0,    0,    0,    0,    0,    0,    0,    0,',
  '20171002, 1200,  , EX, 48.0N,  17.0W,  30, 1016,    0,    0,    0,    0,    0,    0,    0,    0,    0,    0,    0,    0,'])

From above, it is evident that the data is not very structured. We will parse each record into a nicely formatted dictionary for easier analysis.

1b. Format Into Dictionary

hurricanes = {}

for record in records:
    hurricane_id = record[0].split(',')[0]
    hurricane_name = record[0].split(',')[1].strip()
    
    # add to dict
    hurricanes.setdefault(hurricane_id, {'name': hurricane_name})
    hurricanes[hurricane_id].setdefault('datapoints', [])
    
    for datapoint in record[1]:
        data_list = [x.strip() for x in datapoint.split(',')]
        datapoint_dict = {'record_date': data_list[0],
                          'time': data_list[1],
                          'record_identifier': data_list[2], # L = landfall, I = peak intensity
                          'storm_status': data_list[3],
                          'lat': convert_lat_lon(data_list[4], 'lat'), 
                          'lon': convert_lat_lon(data_list[5], 'lon'),
                          'max_wind': data_list[6], 
                          'min_pressure': data_list[7],
                          'category': hurr_category(data_list[6])}
        
        hurricanes[hurricane_id]['datapoints'].append(datapoint_dict)

Lets see what that same record for Hurricane RINA looks like now:

hurricanes['AL192017']

{'name': 'RINA',
 'datapoints': [{'record_date': '20171104',
   'time': '1200',
   'record_identifier': '',
   'storm_status': 'LO',
   'lat': 25.5,
   'lon': -52.3,
   'max_wind': '25',
   'min_pressure': '1013',
   'category': 'TS'},
  {'record_date': '20171104',
   'time': '1800',
   'record_identifier': '',
   'storm_status': 'LO',
   'lat': 26.5,
   'lon': -52.1,
   'max_wind': '25',
   'min_pressure': '1013',
   'category': 'TS'},
  {'record_date': '20171105',
   'time': '0000',
   'record_identifier': '',
   'storm_status': 'LO',
   'lat': 27.5,
   'lon': -52.0,
   'max_wind': '25',
   'min_pressure': '1013',
   'category': 'TS'},
  {'record_date': '20171105',
   'time': '0600',
   'record_identifier': '',
   'storm_status': 'LO',
   'lat': 28.4,
   'lon': -52.0,
   'max_wind': '25',
   'min_pressure': '1013',
   'category': 'TS'},
  {'record_date': '20171105',
   'time': '1200',
   'record_identifier': '',
   'storm_status': 'LO',
   'lat': 29.0,
   'lon': -51.9,
   'max_wind': '25',
   'min_pressure': '1013',
   'category': 'TS'},
  {'record_date': '20171105',
   'time': '1800',
   'record_identifier': '',
   'storm_status': 'TD',
   'lat': 29.2,
   'lon': -51.7,
   'max_wind': '30',
   'min_pressure': '1012',
   'category': 'TS'},
  {'record_date': '20171106',
   'time': '0000',
   'record_identifier': '',
   'storm_status': 'TD',
   'lat': 29.1,
   'lon': -51.2,
   'max_wind': '30',
   'min_pressure': '1011',
   'category': 'TS'},
  {'record_date': '20171106',
   'time': '0600',
   'record_identifier': '',
   'storm_status': 'TD',
   'lat': 29.0,
   'lon': -50.7,
   'max_wind': '30',
   'min_pressure': '1010',
   'category': 'TS'},
  {'record_date': '20171106',
   'time': '1200',
   'record_identifier': '',
   'storm_status': 'TD',
   'lat': 29.1,
   'lon': -50.4,
   'max_wind': '30',
   'min_pressure': '1010',
   'category': 'TS'},
  {'record_date': '20171106',
   'time': '1800',
   'record_identifier': '',
   'storm_status': 'TD',
   'lat': 29.4,
   'lon': -50.2,
   'max_wind': '30',
   'min_pressure': '1009',
   'category': 'TS'},
  {'record_date': '20171107',
   'time': '0000',
   'record_identifier': '',
   'storm_status': 'TS',
   'lat': 30.0,
   'lon': -50.0,
   'max_wind': '35',
   'min_pressure': '1008',
   'category': 'TS'},
  {'record_date': '20171107',
   'time': '0600',
   'record_identifier': '',
   'storm_status': 'TS',
   'lat': 30.8,
   'lon': -49.9,
   'max_wind': '40',
   'min_pressure': '1006',
   'category': 'TS'},
  {'record_date': '20171107',
   'time': '1200',
   'record_identifier': '',
   'storm_status': 'TS',
   'lat': 31.8,
   'lon': -49.5,
   'max_wind': '40',
   'min_pressure': '1004',
   'category': 'TS'},
  {'record_date': '20171107',
   'time': '1800',
   'record_identifier': '',
   'storm_status': 'TS',
   'lat': 33.0,
   'lon': -49.1,
   'max_wind': '40',
   'min_pressure': '1002',
   'category': 'TS'},
  {'record_date': '20171108',
   'time': '0000',
   'record_identifier': '',
   'storm_status': 'TS',
   'lat': 34.6,
   'lon': -48.7,
   'max_wind': '45',
   'min_pressure': '999',
   'category': 'TS'},
  {'record_date': '20171108',
   'time': '0600',
   'record_identifier': '',
   'storm_status': 'TS',
   'lat': 36.4,
   'lon': -48.7,
   'max_wind': '50',
   'min_pressure': '996',
   'category': 'TS'},
  {'record_date': '20171108',
   'time': '1200',
   'record_identifier': '',
   'storm_status': 'TS',
   'lat': 38.3,
   'lon': -48.8,
   'max_wind': '45',
   'min_pressure': '994',
   'category': 'TS'},
  {'record_date': '20171108',
   'time': '1800',
   'record_identifier': '',
   'storm_status': 'TS',
   'lat': 40.1,
   'lon': -49.0,
   'max_wind': '45',
   'min_pressure': '992',
   'category': 'TS'},
  {'record_date': '20171109',
   'time': '0000',
   'record_identifier': '',
   'storm_status': 'TS',
   'lat': 41.8,
   'lon': -48.8,
   'max_wind': '45',
   'min_pressure': '991',
   'category': 'TS'},
  {'record_date': '20171109',
   'time': '0600',
   'record_identifier': '',
   'storm_status': 'LO',
   'lat': 43.6,
   'lon': -48.0,
   'max_wind': '40',
   'min_pressure': '993',
   'category': 'TS'},
  {'record_date': '20171109',
   'time': '1200',
   'record_identifier': '',
   'storm_status': 'LO',
   'lat': 45.5,
   'lon': -47.0,
   'max_wind': '40',
   'min_pressure': '995',
   'category': 'TS'}]}

Much better!

1c. Feature Engineering for Individual Storms

In order to test our hypotheses, we will need some additional data for each hurricane that isn't directly available at the moment, including:

• Indicator if the storm reached hurricane status
• Whether the storm made landfall
• How long the storm lasted

def datapoint_datetime(datapoint):
    """
    For each record, we want a datetime field to calculate
    length of the storm.
    """
    date = datapoint['record_date']
    tm = datapoint['time']
    
    datapoint.setdefault('dt', datetime.datetime.strptime(date+tm, "%Y%m%d%H%M"))
    
    return datapoint

def is_hurricane(storm):
    """
    Did the storm ever make it to hurricane status?
    """
    status_list = []
    for record in storm['datapoints']:
        status_list.append(record['storm_status'])
    
    hurricane = True if 'HU' in status_list else False
    
    return hurricane

def landfall(storm):
    """
    Did the storm make landfall?
    """
    records_on_land = 0
    for record in storm['datapoints']:
        over_land = np.where(record['record_identifier']=='L', 1, 0)
        records_on_land += over_land
    landfall = True if records_on_land > 0 else False
    
    return landfall

def storm_duration(storm):
    """
    How long did the storm last?
    """
    date_list = []
    for record in storm['datapoints']:
        rec_date = record['dt']
        date_list.append(rec_date)

    duration = max(date_list) - min(date_list)
    
    return duration

def get_year(storm):
    """
    Use the first record of the storm to determine the year, 
    since a majority of the analysis is based on how hurricanes
    have changed over the years.
    """
    year = storm['datapoints'][0]['dt'].year
    
    return year

def category_min(storm, cat):
    """
    Did the storm make landfall?
    """
    records_min_cat = 0
    for record in storm['datapoints']:
        if type(record['category'])==int:
            meets_requirement = np.where(record['category']>=cat, 1, 0)
            records_min_cat += meets_requirement
            
    category_min = True if records_min_cat > 0 else False
    
    return category_min

# Apply functions to create features

for key, value in tqdm(hurricanes.items(), total=len(hurricanes.items())):
    for record in value['datapoints']:
        record = datapoint_datetime(record)
    
    value['is_hurricane'] = is_hurricane(value)
    value['landfall'] = landfall(value)
    value['duration'] = storm_duration(value)
    value['year'] = get_year(value)
    value['cat2_or_more'] = category_min(value, 2) 

<ipython-input-20-892d97af01af>:3: TqdmDeprecationWarning: This function will be removed in tqdm==5.0.0
Please use `tqdm.notebook.tqdm` instead of `tqdm.tqdm_notebook`
  for key, value in tqdm(hurricanes.items(), total=len(hurricanes.items())):

Now, let's look Hurricane Rina again to see what we added

rina = hurricanes['AL172017']
print(f"Number of datapoints for Hurricane Rina: {len(rina)}")

Number of datapoints for Hurricane Rina: 7

Great, the data looks ready for some analysis to test our hypotheses! Let's save this cleaned up version before proceeding.

2. Analyze Hurricane Paths

Using our cleaned HURDAT data, we will plot and compare the paths of hurricanes based on different time periods, to get a visual sense of how hurricanes have changed over time due to climate change.\

2a. Compare Hurricane Paths

Before diving into tesitng our hypotheses, we want to look at some maps to visualize how hurricanes have changed over the last 100 years. We will start by looking at the periods 1918-1927 vs 2008-2017.

Our hypothesis states that we expect:

1. More hurricanes in the 2008-2017 period
2. A higher number of very strong storms (Cat 4 & 5)

def get_specific_years(hurricane_data, year_start, year_end):
    """
    Takes an input of the JSON data with hurricane ID as key
    Returns an array of dicts for the range (year_end-1).
    """
    year_range = range(year_start, year_end)

    filtered_data = [storm_data for storm_id, storm_data in hurricane_data.items() if storm_data['year'] in year_range]
    
    return filtered_data

def create_path_map(storm_list, category_colors):
    """
    Given a list of storms to plot and the mapping of category strength
    to color, this function will create a folium plot of all the input
    storms and return the map as a folium object.
    """
    # Initialize map
    hurr_map = folium.Map(location=[30, -70], zoom_start = 4)

    # Plot each storm's path with info in tooltip
    for storm in storm_list:
        storm_name = storm['name']

        for idx in range(len(storm['datapoints'])):
            point = storm['datapoints'][idx]

            # Handle last datapoint which will not have a future datapoint
            if idx+1 < len(storm['datapoints']):
                # Current point
                day_points = [(point['lat'], point['lon']), 
                              (storm['datapoints'][idx+1]['lat'], storm['datapoints'][idx+1]['lon'])]

                hurr_info = f"""
                Storm:        {storm_name}
                Date:         {point['record_date']}
                Time:         {point['time']}
                Status:       {point['storm_status']}
                Category:     {point['category']}
                Max Wind:     {point['max_wind']}mph
                Min Pressure: {point['min_pressure']}millibars
                """
                # Plot line
                folium.PolyLine(day_points, 
                                tooltip=storm_name, 
                                popup=hurr_info, 
                                color=category_colors[point['category']],
                               ).add_to(hurr_map)
                
    return hurr_map

# hurr_1908_1917 = get_specific_years(hurricanes, 1908, 1918)
hurr_1908_1917 = get_specific_years(hurricanes, 1960, 1970)

hurr_2008_2017 = get_specific_years(hurricanes, 2008, 2018)

# Define mapping for colors depending on category
category_colors = {1: '#ffb3b3',
                   2: '#ff8080',
                   3: '#e60000',
                   4: '#ff0000',
                   5: '#580808',
                   'TS': '#ffe6e6'}

map_1908_1917 = create_path_map(hurr_1908_1917, category_colors)
map_2008_2017 = create_path_map(hurr_2008_2017, category_colors)

map_1908_1917

Make this Notebook Trusted to load map: File -> Trust Notebook

map_2008_2017

Make this Notebook Trusted to load map: File -> Trust Notebook

3. Exploratory Data Analysis of Atlantic Hurricanes 1850-2017

After cleaning and parsing the HURDAT dataset on Atlantic Ocean hurricanes, we will do some initial exploratiry data analysis to get a better feel for the data. The aim of this is to potentially uncover some interesting trends and features of the dataset through visualizations and analysis.

3a. Number of Storms

num_hurricanes = len(hurricanes)
print(f'There are data for {num_hurricanes} storms in this dataset.')

There are data for 1848 storms in this dataset.

a = range(1851, 2018)

# Count number of storms
storm_count_by_yr = dict((el,0) for el in a)
for storm_id, storm_data in hurricanes.items():
    year = storm_data['year']
    storm_count_by_yr[year] += 1
    
# Count number of hurricanes
hurricane_count_by_yr = dict((el,0) for el in a)
for storm_id, storm_data in hurricanes.items():
    if storm_data['is_hurricane']==True:
        year = storm_data['year']
        hurricane_count_by_yr[year] += 1

storm_id

'AL192017'

x = list(storm_count_by_yr.keys())
y1 = list(storm_count_by_yr.values())
y2 = list(hurricane_count_by_yr.values())

trace1 = go.Scatter(
    x = x,
    y = y1,
    mode = 'lines+markers',
    name='Number of Tropical Storms'
)

trace2 = go.Scatter(
    x = x,
    y = y2,
    mode = 'lines+markers',
    name='Number of Hurricanes'
)

data = [trace1, trace2]

layout = go.Layout(
    title=go.layout.Title(
        text='Number of Tropical Storms and Hurricanes by Year',
        xref='paper',
        x=0
    ),
    xaxis=go.layout.XAxis(
        title=go.layout.xaxis.Title(
            text='Year',
            font=dict(
                family='Courier New, monospace',
                size=18,
                color='#7f7f7f'
            )
        )
    ),
    yaxis=go.layout.YAxis(
        title=go.layout.yaxis.Title(
            text='Number of Storms',
            font=dict(
                family='Courier New, monospace',
                size=18,
                color='#7f7f7f'
            )
        )
    )
)

fig1 = go.Figure(data=data, layout=layout)
iplot(fig1, filename='storms-by-year')

There seems to be a definite trend towards an increasing number of storms since 1851. Let's go deeper and look at the number of hurricanes in each year. A hurricane is a storm categorized as having wind speed > 74mph.

What's interesting to note is that there are more years with > 10 hurricanes compared to in the past. In the period 1851-1990, there were 9 years where the number of hurricanes was 10 or greater. Since 1990, there have been 6 years with 10 or more hurricanes.

3b. Storm Duration Over Time

According to the HURDAT2 data dictionary, data on tropical storms and depressions were only started in the 1950s and 1960s. If we include all data from 1851, this may skew the results as data from 1850-1950 only includes the major storms

# Get list if durations for each storm system
b = range(1970, 2018)
storm_duration_by_yr = dict((el,[]) for el in b)
for storm_id, storm_data in hurricanes.items():
    year = storm_data['year']
    if year in b:
        duration = storm_data['duration']
        storm_duration_by_yr[year].append(duration)

# using the list, get the average
avg_storm_duration_by_yr = dict((yr,[]) for yr in storm_duration_by_yr.keys())

for year, amounts in storm_duration_by_yr.items():
    day = np.mean(storm_duration_by_yr[year]).days
    day_decimal = np.mean(storm_duration_by_yr[year]).seconds / 86400 # 1 day = 86400 seconds
    avg_time = day + day_decimal
    avg_storm_duration_by_yr[year] = avg_time

x = list(avg_storm_duration_by_yr.keys())
y1 = list(avg_storm_duration_by_yr.values())

trace1 = go.Scatter(
    x = x,
    y = y1,
    mode = 'lines+markers',
    name='Average Duration (Days)'
)

data = [trace1]

layout = go.Layout(
    title=go.layout.Title(
        text='Average Duration of Tropical Storms by Year',
        xref='paper',
        x=0
    ),
    xaxis=go.layout.XAxis(
        title=go.layout.xaxis.Title(
            text='Year',
            font=dict(
                family='Courier New, monospace',
                size=18,
                color='#7f7f7f'
            )
        )
    ),
    yaxis=go.layout.YAxis(
        title=go.layout.yaxis.Title(
            text='Average Duration (Days)',
            font=dict(
                family='Courier New, monospace',
                size=18,
                color='#7f7f7f'
            )
        )
    )
)

fig2 = go.Figure(data=data, layout=layout)
iplot(fig2, filename='avg-duration-by-year')

3c. Number of Landfalls Over Time

# Count number of landfalls
landfall_count_by_yr = dict((el,0) for el in a)
for storm_id, storm_data in hurricanes.items():
    if storm_data['landfall']==True:
        year = storm_data['year']
        landfall_count_by_yr[year] += 1

x = list(landfall_count_by_yr.keys())
y1 = list(landfall_count_by_yr.values())

trace1 = go.Scatter(
    x = x,
    y = y1,
    mode = 'lines+markers',
    name='Number of Landfalls'
)

data = [trace1]

layout = go.Layout(
    title=go.layout.Title(
        text='Number of Tropical Storms and Hurricanes Making Landfall by Year',
        xref='paper',
        x=0
    ),
    xaxis=go.layout.XAxis(
        title=go.layout.xaxis.Title(
            text='Year',
            font=dict(
                family='Courier New, monospace',
                size=18,
                color='#7f7f7f'
            )
        )
    ),
    yaxis=go.layout.YAxis(
        title=go.layout.yaxis.Title(
            text='Number of Landfalls',
            font=dict(
                family='Courier New, monospace',
                size=18,
                color='#7f7f7f'
            )
        )
    )
)

fig3 = go.Figure(data=data, layout=layout)
iplot(fig3, filename='avg-duration-by-year')

4. Hypothesis Analysis

Hypothesis 1

Storms are making more frequent landfalls due to increased strength and duration of systems.

We will focus our analysis on the following metrics to understand hurricane activity across different periods of time.

• Number of storms
• Strong storms (hurricane status)
• Number of landfalls
• Average duration of storms

def get_specific_years(hurricane_data, year_start, year_end):
    """
    Takes an input of the JSON data with hurricane ID as key
    Returns an array of dicts for the range (inclusive).
    """
    year_range = range(year_start, year_end)

    filtered_data = [storm_data for storm_id, storm_data in hurricane_data.items() if storm_data['year'] in year_range]
    
    return filtered_data

def period_metrics(period_name, hurr_list):
    """
    Given a list of hurricanes, this function will calculate metrics for
    that entire list, including number of storms, average duration of storms, 
    number of hurricanes, number of landfall, and values normalized into percentages.
    """
    # 1. Number of storms
    num_storms = len(hurr_list)

    # 2. Average duration of storms
    duration_list = []
    for storm in hurr_list:
        duration_list.append(storm['duration'])
    duration_mean_days = np.mean(duration_list).days

    # 3. Number of hurricanes and storm landfalls
    num_hurricanes = sum([1 for storm in hurr_list if storm['is_hurricane']==True])
    num_landfalls = sum([1 for storm in hurr_list if storm['landfall']==True])
    
    # 4. Normalized percentages for hurricanes landfalls
    pct_landfall = num_landfalls / num_storms
    pct_hurricanes = num_hurricanes / num_storms

    metrics_dict = {'period': period_name,
                'num_storms': num_storms, 
                'avg_duration_days': duration_mean_days, 
                'num_hurricanes': num_hurricanes, 
                'pct_hurricanes': pct_hurricanes, 
                'num_landfalls': num_landfalls, 
                'pct_landfall': pct_landfall}
    
    return metrics_dict

def create_list_of_periods(yr_to_start, yr_to_end, interval):
    """
    Creates an array of dictionaries with the range of the period
    based on the inputted years and interval.
    E.g. for 10 year periods in the years 1917-2017:
        yr_to_start = 1917
        yr_to_end = 2017
        interval = 10
    """
    periods = []
    for yr in list(range(yr_to_start, yr_to_end, interval)):
        start = yr + 1 # to avoid overlapping years across periods
        end = yr + interval
        periods.append({'name': f"{start}-{end}", 'start': start, 'end': end})
        
    return periods

Analyze Historical Averages for 10-year periods

In the 100-year period between 1918-2017, we will analyze the key statistics of each 10-year sub-period (e.g. 1918-1927, 1928-1937, etc).

Obtain history

periods_10yr = create_list_of_periods(1917, 2017, 10)

# Calculate metrics for each period
metrics_10yr = []
for period in periods_10yr:
    hurricane_list = get_specific_years(hurricanes, period['start'], period['end'])
    hurr_latest_metrics = period_metrics(period['name'], hurricane_list)
    metrics_10yr.append(hurr_latest_metrics)

# Convert it into a dataframe to make it easy to read
hurricane_stats_10yr = pd.DataFrame(metrics_10yr)
hurricane_stats_10yr

	period	num_storms	avg_duration_days	num_hurricanes	pct_hurricanes	num_landfalls	pct_landfall
0	1918-1927	63	7	36	0.571429	21	0.333333
1	1928-1937	100	7	48	0.480000	29	0.290000
2	1938-1947	83	5	42	0.506024	29	0.349398
3	1948-1957	120	6	64	0.533333	55	0.458333
4	1958-1967	88	7	53	0.602273	18	0.204545
5	1968-1977	185	5	50	0.270270	1	0.005405
6	1978-1987	142	5	47	0.330986	13	0.091549
7	1988-1997	128	6	55	0.429688	50	0.390625
8	1998-2007	158	7	75	0.474684	76	0.481013
9	2008-2017	140	7	59	0.421429	64	0.457143

hurricane_stats_10yr[['period', 'avg_duration_days']]

	period	avg_duration_days
0	1918-1927	7
1	1928-1937	7
2	1938-1947	5
3	1948-1957	6
4	1958-1967	7
5	1968-1977	5
6	1978-1987	5
7	1988-1997	6
8	1998-2007	7
9	2008-2017	7

Visualize

We will start by visualizing absolute values (number of hurricanes, storms, etc).

sns.set_style("darkgrid", {"axes.facecolor": ".9"})
ax_abs_10yr = sns.lineplot(x="period", y="num_storms", label='Number of Storms', data=hurricane_stats_10yr)
ax_abs_10yr = sns.lineplot(x="period", y="num_hurricanes", label='Number of Hurricanes', data=hurricane_stats_10yr)
ax_abs_10yr = sns.lineplot(x="period", y="num_landfalls", label='Number of Landfalls', data=hurricane_stats_10yr)
ax_abs_10yr.set_title("10-year Period analysis of Hurricanes")
ax_abs_10yr.set(xlabel='Time Period', ylabel='Number of Storms')

for item in ax_abs_10yr.get_xticklabels():
    item.set_rotation(45)

Now, lets visualize the percentage-adjusted metrics (percentage of storms as hurricanes, landfall, etc.) to get a better sense of whether there are more powerful storms now compared to historically.

sns.set_style("darkgrid", {"axes.facecolor": ".9"})
ax_pct_10yr = sns.lineplot(x="period", y="pct_hurricanes", label='% of Storms as Hurricanes', data=hurricane_stats_10yr)
ax_pct_10yr = sns.lineplot(x="period", y="pct_landfall", label='% of Storms Making Landfall', data=hurricane_stats_10yr)
ax_pct_10yr.set_title("10-year Period Analysis of Hurricanes")
ax_pct_10yr.set(xlabel='Time Period', ylabel='Percentage')

for item in ax_pct_10yr.get_xticklabels():
    item.set_rotation(45)

# Look closely at Landfalls
sns.set_style("darkgrid", {"axes.facecolor": ".9"})
ax_land_10yr = sns.lineplot(x="period", y="num_landfalls", label='Number of Landfalls', data=hurricane_stats_10yr, ci=95)
ax_land_10yr.set_title("10-year Period Number of Landfalls")
ax_land_10yr.set(xlabel='Time Period', ylabel='Number of Storms Making Landfalls')

for item in ax_land_10yr.get_xticklabels():
    item.set_rotation(45)

Analyze Historical Averages for 5-year periods

It would be interesting to see if more granular intervals can also capture the general trend of climate change, or whether this range may just be too small to notice anything meaningful due to the variability of hurricanes. This would be 5-year sub-periods (e.g. 1918-1922, 1923-1927, etc).

periods_5yr = create_list_of_periods(1917, 2017, 5)

# Calculate metrics for each period
metrics_list = []
for period in periods_5yr:
    hurricane_list = get_specific_years(hurricanes, period['start'], period['end'])
    hurr_latest_metrics = period_metrics(period['name'], hurricane_list)
    metrics_list.append(hurr_latest_metrics)

hurricane_stats_5yr = pd.DataFrame(metrics_list)
hurricane_stats_5yr

	period	num_storms	avg_duration_days	num_hurricanes	pct_hurricanes	num_landfalls	pct_landfall
0	1918-1922	23	6	15	0.652174	10	0.434783
1	1923-1927	35	7	18	0.514286	11	0.314286
2	1928-1932	27	6	12	0.444444	8	0.296296
3	1933-1937	58	7	30	0.517241	17	0.293103
4	1938-1942	30	6	17	0.566667	12	0.400000
5	1943-1947	42	5	21	0.500000	15	0.357143
6	1948-1952	54	7	32	0.592593	23	0.425926
7	1953-1957	55	6	27	0.490909	29	0.527273
8	1958-1962	45	7	26	0.577778	18	0.400000
9	1963-1967	38	8	24	0.631579	0	0.000000
10	1968-1972	85	5	27	0.317647	1	0.011765
11	1973-1977	81	4	20	0.246914	0	0.000000
12	1978-1982	84	5	26	0.309524	0	0.000000
13	1983-1987	50	5	19	0.380000	13	0.260000
14	1988-1992	62	6	24	0.387097	16	0.258065
15	1993-1997	56	7	27	0.482143	32	0.571429
16	1998-2002	66	6	35	0.530303	30	0.454545
17	2003-2007	78	7	36	0.461538	37	0.474359
18	2008-2012	69	6	30	0.434783	34	0.492754
19	2013-2017	52	7	19	0.365385	23	0.442308

Visualize

We will start by visualizing absolute values (number of hurricanes, storms, etc).

sns.set_style("darkgrid", {"axes.facecolor": ".9"})
ax_abs_5yr = sns.lineplot(x="period", y="num_storms", label='Number of Storms', data=hurricane_stats_5yr)
ax_abs_5yr = sns.lineplot(x="period", y="num_hurricanes", label='Number of Hurricanes', data=hurricane_stats_5yr)
ax_abs_5yr = sns.lineplot(x="period", y="num_landfalls", label='Number of Landfalls', data=hurricane_stats_5yr)
ax_abs_5yr.set_title("5-year Period analysis of Hurricanes")
ax_abs_5yr.set(xlabel='Time Period', ylabel='Number of Storms')

for item in ax_abs_5yr.get_xticklabels():
    item.set_rotation(45)

Now, lets visualize the percentage-adjusted metrics (percentage of storms as hurricanes, landfall, etc.) to get a better sense of whether there are more powerful storms now compared to historically.

sns.set_style("darkgrid", {"axes.facecolor": ".9"})
ax_pct_5yr = sns.lineplot(x="period", y="pct_hurricanes", label='% of Storms as Hurricanes', data=hurricane_stats_5yr)
ax_pct_5yr = sns.lineplot(x="period", y="pct_landfall", label='% of Storms Making Landfall', data=hurricane_stats_5yr)
ax_pct_5yr.set_title("5-year Period Analysis of Hurricanes")
ax_pct_5yr.set(xlabel='Time Period', ylabel='Percentage')

for item in ax_pct_5yr.get_xticklabels():
    item.set_rotation(45)

Look at Average duration of storms

sns.set_style("darkgrid", {"axes.facecolor": ".9"})
ax_duration_10yr = sns.lineplot(x="period", y="avg_duration_days", label='% of Storms Making Landfall', data=hurricane_stats_10yr)
ax_duration_10yr.set_title("10-year Period Analysis of Hurricanes")
ax_duration_10yr.set(xlabel='Time Period', ylabel='Percentage')

for item in ax_duration_10yr.get_xticklabels():
    item.set_rotation(45)

2 year periods

periods_2yr = create_list_of_periods(1917, 2017, 2)

# Calculate metrics for each period
metrics_2yr = []
for period in periods_2yr:
    hurricane_list = get_specific_years(hurricanes, period['start'], period['end'])
    hurr_latest_metrics = period_metrics(period['name'], hurricane_list)
    metrics_2yr.append(hurr_latest_metrics)

# Convert it into a dataframe to make it easy to read
hurricane_stats_2yr = pd.DataFrame(metrics_2yr)
hurricane_stats_2yr

	period	num_storms	avg_duration_days	num_hurricanes	pct_hurricanes	num_landfalls	pct_landfall
0	1918-1919	6	5	4	0.666667	2	0.333333
1	1920-1921	5	6	4	0.800000	3	0.600000
2	1922-1923	5	8	3	0.600000	0	0.000000
3	1924-1925	11	7	5	0.454545	3	0.272727
4	1926-1927	11	9	8	0.727273	3	0.272727
5	1928-1929	6	8	4	0.666667	3	0.500000
6	1930-1931	3	11	2	0.666667	1	0.333333
7	1932-1933	15	7	6	0.400000	4	0.266667
8	1934-1935	13	6	7	0.538462	5	0.384615
9	1936-1937	17	5	7	0.411765	5	0.294118
10	1938-1939	9	5	4	0.444444	3	0.333333
11	1940-1941	9	7	6	0.666667	3	0.333333
12	1942-1943	11	5	4	0.363636	2	0.181818
13	1944-1945	14	5	8	0.571429	4	0.285714
14	1946-1947	7	4	3	0.428571	5	0.714286
15	1948-1949	10	6	6	0.600000	5	0.500000
16	1950-1951	16	8	11	0.687500	6	0.375000
17	1952-1953	11	6	5	0.454545	3	0.272727
18	1954-1955	16	5	7	0.437500	8	0.500000
19	1956-1957	12	5	4	0.333333	5	0.416667
20	1958-1959	12	7	7	0.583333	4	0.333333
21	1960-1961	8	7	4	0.500000	6	0.750000
22	1962-1963	5	8	3	0.600000	0	0.000000
23	1964-1965	12	8	6	0.500000	0	0.000000
24	1966-1967	11	8	7	0.636364	0	0.000000
25	1968-1969	15	5	4	0.266667	0	0.000000
26	1970-1971	19	5	5	0.263158	0	0.000000
27	1972-1973	19	4	3	0.157895	0	0.000000
28	1974-1975	20	4	4	0.200000	0	0.000000
29	1976-1977	21	5	6	0.285714	0	0.000000
30	1978-1979	24	4	5	0.208333	0	0.000000
31	1980-1981	18	6	9	0.500000	0	0.000000
32	1982-1983	8	4	2	0.250000	0	0.000000
33	1984-1985	20	5	5	0.250000	1	0.050000
34	1986-1987	10	5	4	0.400000	2	0.200000
35	1988-1989	19	5	5	0.263158	6	0.315789
36	1990-1991	16	7	8	0.500000	1	0.062500
37	1992-1993	10	6	4	0.400000	2	0.200000
38	1994-1995	12	4	3	0.250000	7	0.583333
39	1996-1997	13	10	9	0.692308	9	0.692308
40	1998-1999	14	9	10	0.714286	7	0.500000
41	2000-2001	19	5	8	0.421053	5	0.263158
42	2002-2003	14	6	4	0.285714	9	0.642857
43	2004-2005	16	8	9	0.562500	8	0.500000
44	2006-2007	10	8	5	0.500000	3	0.300000
45	2008-2009	17	7	8	0.470588	11	0.647059
46	2010-2011	21	7	12	0.571429	11	0.523810
47	2012-2013	19	8	10	0.526316	7	0.368421
48	2014-2015	9	7	6	0.666667	6	0.666667
49	2016-2017	16	8	7	0.437500	9	0.562500

Visualize

We will start by visualizing absolute values (number of hurricanes, storms, etc).

sns.set_style("darkgrid", {"axes.facecolor": ".9"})
ax_abs_2yr = sns.lineplot(x="period", y="num_storms", label='Number of Storms', data=hurricane_stats_2yr)
ax_abs_2yr = sns.lineplot(x="period", y="num_hurricanes", label='Number of Hurricanes', data=hurricane_stats_2yr)
ax_abs_2yr = sns.lineplot(x="period", y="num_landfalls", label='Number of Landfalls', data=hurricane_stats_2yr)
ax_abs_2yr.set_title("10-year Period analysis of Hurricanes")
ax_abs_2yr.set(xlabel='Time Period', ylabel='Number of Storms')
ax_abs_2yr.set_xticklabels(range(1918, 2017))
# ax_abs_2yr.xticks(range(1918, 2017, 10))
for item in ax_abs_2yr.get_xticklabels():
    item.set_rotation(45)  

<ipython-input-70-d16a247c26a8>:7: UserWarning:

FixedFormatter should only be used together with FixedLocator

Hypothesis 2 (TBD)

Storms have increased in strength (wind speed) and intensity (minimum pressure) over time due to rising sea temperatures that have been observed (approximately 0.13 degrees Celsius per decade over the last 100 years).

TODO:

• Anaysis of stats above
• Filter by storm for each state hit